1. Field of the Invention
The present invention relates to a method for determining, from an image of a fractured medium such as, for example, a geologic formation crossed by a network of fractures, hydraulic parameters providing the best large-scale characterization possible thereof, i.e. the large-scale equivalent permeability, the permeability of blocks and the matrix-fracture exchange coefficient .alpha..
2. Description of the Prior Art.
Fractured reservoirs are commonly described as comprising two contrasting mediums, a matrix medium containing the largest part of the oil in place and having a low permeability, and a more or less complex, generally highly permeable medium with fractures characterized by the respective density, length, orientation, inclination and opening thereof.
In a single-phase context, the problem of small-scale well tests simulation is written in the form of a partial differential equation: ##EQU1## where p(r,t) denotes the value of the pressure at point r at the time t, .mu. is the viscosity of the fluid and .phi. the porosity, c.sub.t the total compressibility (rock+fluid) and k(r) the permeability. This equation is to be solved with boundary conditions and initial conditions suited to the N boreholes and to the reservoir boundaries so as to ensure the existence and the uniqueness of the solution thereof, by using as basic data a previously established permeability map k(r). The measurements obtained in the field are the N records of the pressures p(r.sub.i, t), i=1, N, and of the fluid flow rates in these N wells. The practical problem consists in finding information on the permeability map k(r) as a function of these N measurements. All the techniques used to that end are referred to as well test interpretation.
In order to establish such a map in the case of fractured reservoirs, modelling softwares are for example used which generate, from the information available on the medium: photographs, interpreted seismic data, geological surveys, measurements in wells drilled therethrough, core analyses, etc., spatial representations of the network of fractures included in a matrix that is itself permeable according to precise rules described for example in patents French Patents 2,725,795, 2,725,824 and 2,783,073 filed by the applicant. Since the zones of the reservoir around the wells are imperfectly known, these modelling softwares themselves use parameters that are not known well enough. Obtaining of this knowledge from well tests would therefore be desirable. In this fractured context, all the techniques used to that effect are referred to as interpretation of well tests in a fractured medium.
After using this type of modelling software, the oil reservoir .OMEGA. is partitioned into at least two subdomains .OMEGA..sub.f and .OMEGA..sub.m representing the fractures and the matrix respectively. The value of permeability k(r) is k.sub.f (r) if position r is in the fractures subdomain .OMEGA..sub.f, and k.sub.m (r) otherwise, i.e. if r is in matrix .OMEGA..sub.m. Generally, the permeabilities k.sub.f (r) of the fractures are much higher than those, k.sub.m (r), of the matrix: k.sub.f (r)&gt;&gt;k.sub.m (r).
Since the typical dimensions of the fractures (thickness of a few microns or more, very variable lengths) and their great density (several fractures per cubic meter of rock) preclude in practice the use of a fine-grid model for solving partial differential equation (1), it is usual in the petroleum sphere to consider a large-scale model as described for example by:
J. E. Warren & P. J. Root, The Behavior of Natural Fracture Reservoirs, The Society of Petroleum Engineers Journal 3, No. 3, 245-255, 1963. PA1 Quintard M et al, One- and Two-Equation Models for Transient Diffusion Processes in Two-Phase Systems, Advances in Heat Transfer, 23, 369-464, 1993 describes a method allowing obtaining the value of these parameters which comprises solving three "closing problems" defined on a representative portion of the medium. Calculation however requires solving stationary boundary value problems, which requires fine-grid representations of the reservoir. The notion of representative volume is delicate and the practical determination thereof is unknown. Furthermore, the computing cost of these methods becomes large when the network of fractures is complex because they do not spare the necessity of solving large-size linear systems (several million unknowns) if several thousand fractures are to be properly represented. PA1 McCarthy J. F., Continuous-time Random Walks on Random Media, J. Phys. A: Math. Gen. 26, 2495-2503, 1993.
According to this model, which is an idealization of reality, an elementary volume is represented as a set of identical parallelepipedic blocks limited by a network of continuous uniform orthogonal fractures oriented in the direction of one of the three principal axes. The fluid flow through the reservoir occurs through the fracture medium only and fluid exchanges appear locally between the fractures and the matrix blocks. The following equations are solved: ##EQU2##
Pressures P.sub.m and P.sub.f are large-scale mean values of the pressures in the matrix and in the fracture. Coefficients K.sub.m and K.sub.f represent the large-scale equivalent permeabilities of the matrix and of the network of fractures. Quantities .phi..sub.f and .phi..sub.m denote the volume proportions of fractures and of the matrix. We have of course the following relation: .phi..sub.f +.phi..sub.m =1. Term q.sub.m (r, t) represents a source term of exchange between the two pressures which, in the "pseudosteady" case, is written as follows: ##EQU3##
This term represents the mass exchange between the two mediums, which justifies the name "term of exchange" and the name "exchange coefficient" applied to .alpha.. It can be noted that quantity .alpha.K.sub.m /.phi..mu.c.sub.t is homogeneous to the inverse of a time.
It is essential to understand that, at the scale considered, individual fractures have disappeared from the description. In particular, parameters K.sub.m, K.sub.f and .alpha. are generally independent of the point r considered, and details of the properties of the network of fractures appear only via the value of the parameters.
The physical idea at the root of the present modelling is to consider that the fluid stored in the matrix moves towards the wells via the network of fractures. In the case of a medium consisting of cubic matrix blocks of side L separated by fractures, .alpha. is of the order of 1/L.sup.2, whence the expression "equivalent block size" that is sometimes encountered to designate this coefficient. Other more complex representations of term q.sub.m (r, t) exist in the form of temporal convolutions and are referred to as "transient models".
This formulation, notably known from Warren et al, mentioned above, when it is extended to multiphase flows, represents what is referred to as a "double medium" description. It serves as a basis for all the simulations of flows in fractured reservoirs that are usually performed in the petroleum sphere. In the particular case of single-phase well tests, knowledge of the analytical solutions to Eqs.(2) allows to adjust the parameters from the pressure variations observed in the well(s). In particular, if the well has a radius r.sub.w, coefficient .lambda.=.alpha.K.sub.m (L/r.sub.w).sup.2 is conventionally introduced. In practice, correct adjustments are obtained for well tests, but if it is desired to adjust the description to the field scale, other values can be found for the coefficients.
This model, whose great quality lies in the simplicity thereof, has a great deficiency: parameters K.sub.m, K.sub.f and .alpha. are a priori phenomenological and adjustable. Their relation with the detailed structure of the medium and the flow process is not widely known and few efficient computing tools are available.
Calculation of these parameters requires scaling of the detailed description of the medium. For example:
A technique referred to as "random walk" technique is used to determine the permeability variations of the homogeneous and heterogeneous medium. It is notably described by:
The principle of the random walk technique follows: a certain number of independent fluid particles initially placed on certain grids of a regular grid pattern discretizing the medium studied is considered, and the displacement thereof on the neighbouring grids is studied with a weighted equal probability in order to take account of the medium heterogeneities. Return to a global vision of the medium is performed by calculating the average of the squares of the displacements of all these particles. It is supposed that the problem of the solution of Eq.(1) has been previously discretized by means of a finite-difference or a finite-element method with a regular spatial grid pattern of N grids. The discretized equations are written as follows: ##EQU4## where mj denotes the grids next to grid i.
The N quantities P.sub.i (t) thus represent the pressure evaluated at the center of the N nodes. Coefficients T.sub.ij are the transmissivities between grids i and j, and they are given by the permeability harmonic means k.sub.i and k.sub.j of the two grids considered. The heterogeneous or fractured nature of the medium is thus taken into account via the values of these transmissivities. Coefficient Vi=(.DELTA.X).sup.D is the volume of the grids which is assumed not to depend on i (regular grid pattern hypothesis).
The random walk is considered when as follows. At the initial time t=0 ("zero" iteration, k=0), a grid with subscript i.sub.0 is taken at random and it is assumed that a particle is released therein. At the iteration k of the algorithm, it is assumed that the particle has reached the grid marked by subscript i.sub.k, and that the time is t.sub.k. In order to find its position i.sub.k +1 at iteration k+1, a grid next to i.sub.k and bearing subscript j is selected at random with the following probability: ##EQU5##
The sum on subscript j concerns the same neighbouring grids as in Eq.(4). The probabilities b are normalized, the value of their sum being 1 by construction. The particle has to jump from one site to the other without it being possible for it to remain on the spot. In the particular case where the particle is on a site at the border of the domain, the periodic boundary conditions are used: everything goes off as if the particle could "go out" of the medium and return thereto on the other side (the absolute position thereof being however kept).
This selection being made, the position i.sub.k+1 at iteration k+1 is known, and the timer is then incremented as follows: ##EQU6## where ran(k+1) is a random number of uniform law on interval (0, 1). The logarithm still being negative, the timer increases with subscript k. .DELTA.X is the spatial pitch of the grid pattern. The quantity factorized is homogeneous to a time and represents a typical diffusion time over the distance .DELTA.X. The increase of the timer corresponds to the distribution of the residence times of a particle on the grid considered.
The process can now be iterated and it is stopped when the time t.sub.k exceeds a set value T.sub.max set by the user.
A random path or walk thus consists of a table (of displacements) of size K(T.sub.max) (size depending on the realization) comprising the values of K(T.sub.max) positions i.sub.k and of K(T.sub.max) associated times t.sub.k.
Consider now a great number of such independent random walks (or paths), therefore Q lists of positions and times. The probability for a particle of being on grid i between t and t+dt defines a probability density P.sub.i (t) which could be calculated by means of a very great number of walks.
It can be demonstrated that all of the densities P.sub.i (t) follow Eqs. (5). In the case considered, where the particles start at t=0 from a random site of the network, the initial condition is P.sub.i (t=0)=1/N. If it is decided to always start from the same site i.sub.0, at t=0 we would have P.sub.i (t=0)=.delta.i i.sub.0, where .delta. is the Kroenecker symbol. The discrete Green function of the problem is obtained, i.e. the probability of being on grid "i" at the time t, knowing that we were on grid i.sub.0 at the time t=0.
This type of displacement is referred to as "continuous-time random walk" (CTRW). It is the apparently complex management of the time interval (Eq.(6)) which allows the particles to move at each interval without ever pausing on a grid and precision obtaining of good numerical performances (since each pause would correspond to numerical iterations and therefore to unnecessarily spent computer time).